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Complete  Configuration  Aero-Structural 
Optimization  Using  a  Coupled  Sensitivity 

Analysis  Method 

Joaquim  R.  R.  A.  Martins*  Juan  J.  Alonso  t 
Stanford  University,  Stanford,  CA  94305 

James  J.  Reuther  + 

NASA  Ames  Research  Center,  Moffett  Field,  CA  94035 

This  paper  focuses  on  the  demonstration  of  a  new  integrated  aero-structural  design 
method  for  aerospace  vehicles.  The  approach  combines  an  aero-structural  analysis  solver, 
a  coupled  aero-structural  adjoint  solver,  a  geometry-based  analysis  and  design  integration 
strategy,  and  an  efficient  gradient-based  optimization  algorithm.  The  aero-structural 
solver  ensures  highly  accurate  solutions  by  using  high-fidelity  models  for  both  disciplines 
as  well  as  a  high-fidelity  coupling  procedure.  The  Euler  equations  are  solved  for  the 
aerodynamics  and  a  detailed  finite  element  model  is  used  for  the  primary  structure.  The 
coupled  aero-structural  adjoint  solution  is  used  to  calculate  the  needed  sensitivities  of 
aerodynamic  and  structural  cost  functions  with  respect  to  both  aerodynamic  shape  and 
structural  variables.  The  geometric  outer  mold  line  (OML)  serves  not  only  as  an  interface 
between  the  two  disciplines  for  both  the  state  and  costate  systems,  but  also  as  an  interface 
between  the  numerical  optimization  algorithm  and  the  high-fidelity  analyses.  Another  set 
of  design  variables  parameterizes  a  structure  of  fixed  topology.  Kreisselmeier— Steinhauser 
functions  are  used  to  reduce  the  number  of  structural  constraints  in  the  problem.  Sample 
results  comparing  a  fully  coupled  aero-structural  design  with  a  more  traditional  sequential 
optimization  are  presented. 


Introduction 

During  the  past  decade  the  advancement  of  numer¬ 
ical  methods  for  the  analysis  of  complex  engineer¬ 
ing  problems  such  as  those  found  in  fluid  dynamics 
and  structural  mechanics  has  reached  a  mature  stage: 
many  difficult  numerically  intensive  problems  are  now 
readily  solved  with  modern  computer  facilities.  In  fact, 
the  aircraft  design  community  is  increasingly  using 
computational  fluid  dynamics  (CFD)  and  computa¬ 
tional  structural  mechanics  (CSM)  tools  to  replace 
traditional  approaches  based  on  simplified  theories  and 
wind  tunnel  testing.  With  the  advancement  of  these 
numerical  analysis  methods  well  underway,  the  focus 
for  engineers  is  shifting  toward  integrating  these  anal¬ 
ysis  tools  into  numerical  design  procedures. 

Despite  revolutionary  accomplishments  in  single¬ 
discipline  applications  [23,  24],  progress  towards  the 
development  of  high-fidelity,  multidisciplinary  design 
optimization  (MDO)  methods  has  been  slow.  The 
level  of  coupling  between  disciplines  is  highly  problem 
dependent  and  significantly  affects  the  choice  of  algo¬ 
rithm.  Multiple  difficulties  also  arise  from  the  wide 
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variety  of  design  problems:  an  approach  that  is  ap¬ 
plicable  to  one  problem  may  not  be  compatible  with 
another. 

An  important  feature  that  characterizes  the  various 
solution  strategies  for  MDO  problems  is  the  allowable 
level  of  disciplinary  autonomy  in  the  analysis  and  op¬ 
timization  components.  Excellent  discussions  of  these 
issues  are  presented  by  Sobieski  and  Haftka  [28],  and 
Alexandrov  and  Lewis  [4] .  The  allowable  level  of  disci¬ 
plinary  autonomy  is  usually  inversely  proportional  to 
the  bandwidth  of  the  interdisciplinary  coupling.  Thus, 
for  highly  coupled  problems  it  may  be  necessary  to  re¬ 
sort  to  fully  integrated  MDO,  while  for  more  weakly 
coupled  problems,  modular  strategies  may  hold  an  ad¬ 
vantage  in  terms  of  ease  of  implementation.  With 
these  constraints  in  mind,  a  number  of  ideas  for  solving 
complex  MDO  problems  have  been  developed.  These 
ideas  include  multilevel  optimization  strategies  [3,  15], 
collaborative  optimization  [6,  17,  9],  individual  disci¬ 
pline  feasible  methods  [8],  as  well  as  tightly  coupled 
optimization  procedures.  The  main  difference  between 
the  different  MDO  strategies  is  the  degree  of  coupling 
that  is  required  between  the  disciplines  in  both  the 
analysis  and  the  optimization  procedures. 

In  the  particular  case  of  high-fidelity  aero-structural 
optimization,  the  coupling  between  disciplines  has  a 
very  high  bandwidth.  Furthermore,  the  values  of  the 
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objective  functions  and  constraints  depend  on  highly 
coupled  multidisciplinary  analyses  (MDA).  As  a  result, 
we  believe  that  a  tightly  coupled  MDO  environment  is 
more  appropriate  for  aero-structural  optimization. 

Another  important  consideration  when  selecting 
a  numerical  design  strategy  is  to  choose  between 
gradient-  and  non-gradient-based  approaches.  Some 
combination  of  these  two  approaches  is  likely  to  emerge 
as  the  leading  alternative  in  which  discontinuities  in 
the  design  space,  such  as  changes  in  the  structural 
topology,  are  treated  with  a  non-gradient-based  proce¬ 
dure,  while  efficient  convergence  to  a  local  minimum  is 
achieved  with  a  gradient  approach.  The  present  paper 
employs  a  gradient-based  strategy  that  enables  the  use 
of  hundreds  or  even  thousands  of  design  parameters 
to  achieve  near  optimal  shape-structure  combinations. 
This  level  of  design  detail  —  which  is  arguably  nec¬ 
essary  at  the  end  of  the  preliminary  design  stage  ■ 
cannot  be  treated  by  non-gradient  methods  when  the 
analyses  involve  the  solution  of  a  high-fidelity  aero- 
structural  system. 

In  contrast  with  emerging  single-discipline  design 
methodologies  such  as  aerodynamic  shape  optimiza¬ 
tion  or  structural  optimization,  aero-structural  design 
has  traditionally  been  carried  out  in  a  cut-and-try  ba¬ 
sis.  Aircraft  designers  have  a  pre-conceived  idea  of 
the  shape  of  an  “optimal”  load  distribution  and  then 
tailor  the  jig  shape  of  the  structure  so  that  the  de¬ 
flected  wing  shape  under  a  1-g  load  gives  the  desired 
load  distribution.  While  this  approach  may  suffice  for 
conventional  transport  aircraft,  where  there  is  consid¬ 
erable  accumulated  experience,  in  the  case  of  either 
new  planform  concepts  or  new  flight  regimes,  the  lack 
of  experience  combined  with  the  complexities  of  aero- 
structural  interactions  can  lead  to  designs  that  are  far 
from  optimal. 

This  is  certainly  the  case  in  the  design  of  both  small 
and  large  supersonic  transports,  where  simple  beam 
theory  models  of  the  wing  cannot  be  used  to  accu¬ 
rately  describe  the  behavior  of  the  wing  structure.  In 
some  cases,  these  aircraft  must  even  cruise  for  signifi¬ 
cant  portions  of  their  flight  at  different  Mach  numbers. 
In  addition,  a  variety  of  studies  show  that  supersonic 
transports  exhibit  a  range  of  undesirable  aeroelastic 
phenomena  due  to  the  low  bending  and  torsional  stiff¬ 
ness  that  result  from  wings  with  low  thickness  to  chord 
ratio.  These  phenomena  can  only  be  suppressed  when 
aero-structural  interactions  are  taken  into  account  at 
the  preliminary  design  stage  [5]. 

Unfortunately,  the  modeling  of  the  participating  dis¬ 
ciplines  in  most  of  the  work  that  has  appeared  so  far 
has  remained  at  a  relatively  low  level.  While  use¬ 
ful  at  the  conceptual  design  stage,  lower-order  mod¬ 
els  cannot  accurately  represent  a  variety  of  nonlin¬ 
ear  phenomena  such  as  wave  drag,  which  can  play 
an  important  role  in  the  search  for  the  optimum  de¬ 
sign.  An  exception  to  low- fidelity  modeling  is  the 


recent  work  by  Giunta  [11]  and  by  Maute  et  al.  [22] 
where  aero-structural  sensitivities  are  calculated  using 
higher-fidelity  models. 

The  ultimate  objective  of  our  work  is  to  develop 
an  MDO  framework  for  high-fidelity  analysis  and  op¬ 
timization  of  aircraft  configurations.  The  framework 
is  built  upon  prior  work  by  the  authors  on  aero- 
structural  high-fidelity  sensitivity  analysis  [25,  18,  19]. 
The  objective  of  this  paper  is  to  present  the  current 
capability  of  this  framework  and  to  demonstrate  it  by 
performing  a  simplified  aero-structural  design  of  a  su¬ 
personic  business  jet  configuration. 

This  paper  presents  a  tightly  coupled  approach  to 
high-fidelity  aero-structural  MDO  that  uses  CFD  and 
CSM.  The  following  sections  begin  with  the  descrip¬ 
tion  of  the  aircraft  optimization  problem  we  propose 
to  solve.  We  then  introduce  the  general  formulation  of 
the  sensitivity  equations  followed  by  the  specific  case 
of  the  adjoint  equations  for  the  aero-structural  sys¬ 
tem.  A  detailed  study  of  the  accuracy  and  efficiency 
of  the  aeroelastic  sensitivity  information  is  also  pre¬ 
sented  for  validation  purposes.  Finally,  we  present 
results  of  the  application  of  our  sensitivity  analysis 
method  to  the  full  aero-structural  optimization  of  a 
supersonic  business  jet  and  compare  the  results  with 
the  more  traditional  approach  of  sequential  discipline 
optimizations  where  we  highlight  the  fact  that  only 
closely  coupled  optimization  frameworks  can  arrive  at 
the  true  optimum  of  the  system. 


Aircraft  Optimization  Problem 


For  maximum  lift-to-drag  ratio,  it  is  a  well-known 
result  from  classical  aerodynamics  that  a  wing  must 
exhibit  an  elliptic  lift  distribution.  For  aircraft  design, 
however,  it  is  usually  not  the  lift-to-drag  ratio  we  want 
to  maximize  but  an  objective  function  that  reflects  the 
overall  mission  of  the  particular  aircraft.  Consider,  for 
example,  the  Breguet  range  formula  for  jet-powered 
aircraft 


Range 


V  CL  ^  Wi 
c  CD  ln  Wf  * 


(1) 


where  V  is  the  cruise  velocity  and  c  is  the  thrust  spe¬ 
cific  fuel  consumption  of  the  powerplant.  Cl/Cd  is 
the  ratio  of  lift  to  drag,  and  Wi/Wf  is  the  ratio  of 
initial  and  final  cruise  weights  of  the  aircraft. 

The  Breguet  range  equation  expresses  a  trade-off  be¬ 
tween  the  drag  and  the  empty  weight  of  the  aircraft 
and  constitutes  a  reasonable  objective  function  to  use 
in  aircraft  design.  If  we  were  to  parameterize  a  de¬ 
sign  with  both  aerodynamic  and  structural  variables 
and  then  maximize  the  range  for  a  fixed  initial  cruise 
weight,  subject  to  stress  constraints,  we  would  obtain 
a  lift  distribution  similar  to  the  one  shown  in  Figure  1. 


This  optimum  lift  distribution  trades  off  the  drag 
penalty  associated  with  unloading  the  tip  of  the  wing, 
where  the  loading  contributes  most  to  the  maximum 
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tion. 


Fig.  2:  Natural  laminar  flow  supersonic  business  jet  con¬ 
figuration. 


stress  at  the  root  of  the  wing  structure,  in  order  to 
reduce  the  weight.  The  end  result  is  an  increase  in 
range  when  compared  to  the  elliptically  loaded  wing 
that  results  from  an  increased  weight  fraction  Wt/Wf. 
The  result  shown  in  Figure  1  illustrates  the  need  for 
taking  into  account  the  coupling  of  aerodynamics  and 
structures  when  performing  aircraft  design. 

The  aircraft  configuration  used  in  this  work  is  the 
supersonic  business  jet  shown  in  Figure  2.  This  con¬ 
figuration  is  being  developed  by  the  ASSET  Research 
Corporation  and  is  designed  to  achieve  a  large  percent¬ 
age  of  laminar  flow  on  the  low-sweep  wing,  resulting 
in  decreased  friction  drag  [16].  The  aircraft  is  to  fly  at 
Mach  1.5  and  have  a  range  of  5,300  nautical  miles. 

Detailed  mission  analysis  for  this  aircraft  has  deter¬ 
mined  that  one  count  of  drag  (A Cd  =  0.0001)  is  worth 
310  pounds  of  empty  weight.  This  means  that  to  op¬ 
timize  the  range  of  the  configuration  we  can  minimize 
the  objective  function 

I  =  aCD  +  0W,  (2) 

where  Cd  is  the  drag  coefficient,  W  is  the  structural 
weight  in  pounds  and  a/0  =  3.1  x  106. 

The  aircraft  design  is  parameterized  using  two  types 
of  design  variables.  The  first  type  of  variable  modifies 
the  OML  of  the  configuration  while  the  second  type  of 
variable  controls  the  sizing  of  underlying  structure. 


To  perform  gradient-based  optimization,  we  need 
the  sensitivities  of  the  objective  function  (2)  with  re¬ 
spect  to  all  the  design  variables.  Since  this  objective 
function  is  a  linear  combination  of  the  drag  coefficient 
and  the  structural  weight,  its  sensitivity  can  be  written 
as 


dJ 

dxn 


=  a 


d Cd 

d  Xn 


■P 


d  W 

d  Xn 


(3) 


The  sensitivity  of  the  structural  weight,  dW  /  dxn ,  is 
trivial,  since  the  weight  calculation  is  independent  of 
the  aero-structural  solution.  This  gradient  is  calcu¬ 
lated  analytically  for  the  structural  thickness  variables 
and  by  finite  differences  for  the  OML  variables.  The 
drag  coefficient  sensitivity,  dCo/ dxn,  is  not  trivial 
since  it  does  depend  on  the  aero-structural  solution. 
Details  of  the  methodology  used  to  compute  coupled 
sensitivities  are  presented  in  the  next  section. 

The  OML  design  variables  can  be  applied  to  any 
of  the  components  used  to  define  the  aircraft  geome¬ 
try.  For  each  wing-like  component  (main  wing,  canard, 
horizontal  tail,  etc.),  the  shape  is  modified  at  a  number 
of  pre-specified  airfoil  sections.  Each  of  these  sections 
may  be  independently  modified  while  the  spanwise 
resolution  can  be  controlled  by  the  number  and  po¬ 
sition  of  the  sections.  The  shape  modifications  to 
the  airfoils  are  linearly  lofted  between  stations.  Var¬ 
ious  types  of  design  variables  may  be  applied  to  the 
airfoils:  twist,  leading  and  trailing  edge  droop,  and 
Hicks-Henne  bump  functions,  among  others.  The 
Hicks-Henne  functions  are  of  the  form 


HO  =  Xn 


sin 


t2 


(4) 


where  ti  is  the  location  of  the  maximum  of  the  bump  in 
the  range  0<£<lat£  =  ii,  since  the  maximum  oc¬ 
curs  when  £“  =  1/2,  where  a  =  log(l/2)/logti.  The 
parameter  t2  controls  the  width  of  the  bump.  The 
advantage  of  these  functions  is  that  when  they  are  ap¬ 
plied  to  a  smooth  airfoil,  that  airfoil  remains  smooth. 

The  structural  design  variables  are  the  thicknesses 
of  the  structural  finite  elements.  The  topology  of  the 
structure  remains  unchanged,  i.e.,  the  number  of  spars 
and  ribs  and  their  position  are  fixed  throughout  the 
optimization.  However,  because  the  OML  determines 
the  location  of  the  nodes  of  the  structural  model,  vari¬ 
ations  of  the  OML  have  an  effect  on  the  depth  of  the 
spars  and  ribs  of  the  wing  box. 

Among  the  constraints  to  be  imposed,  the  most  ob¬ 
vious  one  is  that  during  cruise  the  lift  must  always 
equal  the  weight  of  the  aircraft.  In  our  optimization 
problem  we  constrain  the  Cl  by  periodically  adjusting 
the  angle-of-attack  of  the  aircraft  within  the  aero- 
structural  solution  until  the  desired  lift  is  obtained 
within  a  pre-specified  tolerance.  Otherwise,  OML  de¬ 
sign  changes  would  quickly  result  in  lower  drag  coef¬ 
ficients  simply  because  of  reduced  lift.  Some  design 
problems  may  require  that  the  objective  function  be 
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minimized  over  a  range  of  operating  conditions  (mul¬ 
tipoint  design).  In  these  situations,  the  appropriate 
lift  constraint  would  be  imposed  at  each  design  point 
using  the  same  procedure. 

In  addition  to  maintaining  the  Cl,  the  stresses  are 
also  constrained  so  that  the  yield  stress  of  the  mate¬ 
rial  is  not  exceeded  at  a  number  of  load  conditions. 
There  are  typically  thousands  of  finite  elements  de¬ 
scribing  the  structure  of  the  aircraft,  and  it  can  become 
computationally  very  costly  to  treat  these  constraints 
separately.  The  difficulty  of  the  problem  is  that  even 
though  there  are  efficient  ways  of  computing  sensitiv¬ 
ities  of  a  few  functions  with  respect  to  many  design 
variables,  and  of  computing  sensitivities  of  many  func¬ 
tions  with  respect  to  a  few  design  variables,  there  is 
no  known  efficient  method  for  computing  sensitivities 
of  many  functions  with  respect  to  many  design  vari¬ 
ables.  Thus,  we  are  left  to  choose  between  treating 
a  large  number  of  design  variables  or  being  able  to 
handle  multiple  cost  functions  and  constraints. 

In  the  case  of  most  structural  design  problems,  the 
preferred  approach  is  the  direct  method,  where  one  cal¬ 
culates  the  sensitivities  of  the  stress  in  each  element 
independently,  while  limiting  the  number  of  design 
parameters  to  lower  the  total  cost  of  computing  sen¬ 
sitivities.  The  number  of  sensitivity  analyses  needed 
to  compute  the  complete  sensitivity  matrix  scales  in¬ 
dependently  of  the  number  of  stress  constraints  but 
linearly  with  respect  to  the  number  of  design  param¬ 
eters.  Unfortunately,  the  most  efficient  approach  for 
aerodynamics  problems  is  quite  different:  the  number 
of  aerodynamic  cost/constraint  functions  is  relatively 
small,  while  the  number  of  necessary  shape  design 
parameters  is  very  large.  It  is  then  more  efficient 
to  compute  gradients  via  the  adjoint  method,  whose 
cost  scales  linearly  with  the  number  of  cost/constraint 
functions  but  is  independent  of  the  number  of  design 
parameters.  The  coupled  aero-structural  problem  re¬ 
quires  a  compromise  between  these  two  approaches. 
However,  without  a  viable  strategy  to  reduce  the  need 
for  a  large  number  of  aerodynamic  design  variables  we 
are  left  with  the  option  of  trying  to  reduce  the  number 
of  independent  structural  constraints. 

For  this  reason,  we  lump  the  individual  element 
stresses  using  Kreisselmeier - -Steinhauser  (KS)  func¬ 
tions.  In  the  limit,  all  element  stress  constraints  can 
be  lumped  into  a  single  KS  function,  thus  minimizing 
the  cost  of  a  large-scale  aero-structural  design  cycle. 
Suppose  that  we  have  the  following  constraint  for  each 
structural  finite  element, 

9m  =  1  -  —  >  0,  (5) 

ay 

where  am  is  the  element  von  Mises  stress  and  ay  is 
the  yield  stress  of  the  material.  The  corresponding  KS 


function  is  defined  as 


This  function  represents  a  lower  bound  envelope  of  all 
the  constraint  inequalities  and  p  is  a  positive  parame¬ 
ter  that  expresses  how  close  this  bound  is  to  the  actual 
minimum  of  the  constraints.  This  constraint  lumping 
method  is  conservative  and  may  not  achieve  the  same 
optimum  that  a  problem  treating  the  constraints  sep¬ 
arately  would.  However,  the  use  of  KS  functions  has 
been  demonstrated  and  it  constitutes  a  viable  alter¬ 
native,  being  effective  in  optimization  problems  with 
thousands  of  constraints  [2] . 

Having  defined  our  objective  function,  design  vari¬ 
ables  and  constraints,  we  can  now  summarize  the  air¬ 
craft  design  optimization  problem  as  follows: 

minimize  /  =  olCb  +  f3W 
xn  €  M" 

subject  to  Cl  =  Clt 
KS  >  0 

Xn  >  Xnm in- 

The  stress  constraints  in  the  form  of  KS  functions  must 
be  enforced  by  the  optimizer  for  aerodynamic  loads 
corresponding  to  a  number  of  flight  and  dynamic  load 
conditions.  Finally,  a  minimum  gage  is  specified  for 
each  structural  element  thickness. 

Analytic  Sensitivity  Analysis 

Our  main  objective  is  to  calculate  the  sensitivity  of 
a  multidisciplinary  function  of  interest  with  respect  to 
a  number  of  design  variables.  The  function  of  inter¬ 
est  can  be  either  the  objective  function  or  any  of  the 
constraints  specified  in  the  optimization  problem.  In 
general,  such  functions  depend  not  only  on  the  design 
variables,  but  also  on  the  physical  state  of  the  multi¬ 
disciplinary  problem.  Thus  we  can  write  the  function 
as 

I  I(xn,  Vi)i  (7) 

where  x„  represents  the  vector  of  design  variables  and 
Hi  is  the  state  variable  vector. 

For  a  given  vector  xn,  the  solution  of  the  governing 
equations  of  the  multidisciplinary  system  yields  a  vec¬ 
tor  i u,  thus  establishing  the  dependence  of  the  state  of 
the  system  on  the  design  variables.  We  denote  these 
governing  equations  by 

'R-k  (xn,  yi  ( xn ))  —  0.  (8) 

The  first  instance  of  xn  in  the  above  equation  indicates 
the  fact  that  the  residual  of  the  governing  equations 
may  depend  explicitly  on  xn.  In  the  case  of  a  CFD 
solver,  for  example,  changing  the  surface  shape  results 
in  different  values  of  the  residual  for  at  least  the  mesh 
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points  closest  to  the  surface,  even  if  the  solution  is  not 
recomputed.  By  solving  the  governing  equations  we 
determine  the  state,  yi ,  which  depends  implicitly  on 
the  design  variables  through  the  solution  of  the  system. 
These  equations  may  be  nonlinear,  in  which  case  the 
usual  procedure  is  to  drive  residuals,  TZk,  to  zero  using 
an  iterative  method. 

Since  the  number  of  equations  must  equal  the  num¬ 
ber  of  state  variables,  the  ranges  of  the  indices  i  and 
k  are  the  same,  i.e. ,  i,k  =  1, . . . ,  Nn.  In  the  case 
of  a  structural  solver,  for  example,  Nn  is  the  num¬ 
ber  of  unconstrained  degrees  of  freedom,  while  for  a 
CFD  solver,  N-ji  is  the  number  of  mesh  points  multi¬ 
plied  by  the  number  of  state  variables  stored  at  each 
point.  In  the  more  general  case  of  a  multidisciplinary 
system,  TZk  represents  all  the  governing  equations  of 
the  different  disciplines,  including  their  coupling. 


Fig.  3:  Schematic  representation  of  the  governing  equa¬ 
tions  (IZk  =  0),  design  variables  ( x„ ),  state  variables  (yi), 
and  objective  function  (I),  for  an  arbitrary  system. 


A  graphical  representation  of  the  system  of  govern¬ 
ing  equations  is  shown  in  Figure  3,  with  the  design 
variables  xn  as  the  inputs  and  I  as  the  output.  The 
two  arrows  leading  to  /  illustrate  the  fact  that  the 
objective  function  typically  depends  on  the  state  vari¬ 
ables  and  may  also  be  an  explicit  function  of  the  design 
variables. 

As  a  first  step  toward  obtaining  the  derivatives  that 
we  ultimately  want  to  compute,  we  use  the  chain  rule 
to  write  the  total  sensitivity  of  I  as 

d/  dl  dl  d? /i 

—  =  w - fw — -t — )  (9) 

dxn  (JXn  oyi  o.xn 


for  i  =  1, . . . ,  Nn,  n  =  1, . . . ,  Nx.  Index  notation  is 
used  to  denote  the  vector  dot  products.  It  is  impor¬ 
tant  to  distinguish  the  total  and  partial  derivatives  in 
this  equation.  The  partial  derivatives  can  be  directly 
evaluated  by  simply  varying  the  denominator  and  re¬ 
evaluating  the  function  in  the  numerator  while  keep¬ 
ing  everything  else  constant.  The  total  derivatives, 
however,  require  the  solution  of  the  multidisciplinary 
problem.  Thus,  all  the  terms  in  the  total  sensitiv¬ 
ity  equation  (9)  are  inexpensively  computed  except  for 
d  yi/  dxn. 

Since  the  governing  equations  must  always  be  sat¬ 
isfied,  the  total  derivative  of  the  residuals  (8)  with 
respect  to  any  design  variable  must  also  be  zero.  Ex¬ 
panding  the  total  derivative  of  the  governing  equations 
with  respect  to  the  design  variables  we  can  write, 

d TZk  _  dTZk  |  dTZk  d yi  _  Q 

dxn  dxn  dyt  dxn 


for  all  i,k  =  1, . . . ,  N-ji  and  n  =  1, . . . ,  Nx .  This  ex¬ 
pression  provides  us  with  the  means  to  compute  the 
total  sensitivity  of  the  state  variables  with  respect  to 
the  design  variables.  By  rewriting  equation  (10)  as 

dTZk  d yj  _  dTZk  ^ 

<>y,  d.r„  dxn  ’ 

we  can  solve  for  dpi/ dxn  and  substitute  this  result 
into  the  total  derivative  equation  (9),  to  obtain 


dl 

dxn 


—  d yi  /  dxn 

/ - -*■ - 

d I_ 

dyi 

-'S'k 


dTZk  dTZk 
dyi  dxn 


(12) 


The  inversion  of  the  Jacobian  dTZk /dpi  is  not  necessar¬ 
ily  calculated  explicitly.  In  the  case  of  large  iterative 
problems  neither  this  matrix  nor  its  factorization  are 
usually  stored  due  to  their  prohibitive  size. 

The  approach  where  we  first  calculate  dyi  /  dxn  us¬ 
ing  equation  (11)  and  then  use  the  result  in  the  expres¬ 
sion  for  the  total  sensitivity  (12)  is  called  the  direct 
method.  Note  that  solving  for  dyi/  dxn  requires  the 
solution  of  the  matrix  equation  (11)  for  each  design 
variable  xn.  A  change  in  the  design  variable  affects 
only  the  right-hand  side  of  the  equation,  so  for  prob¬ 
lems  where  the  matrix  dTZk/dpi  can  be  explicitly  fac¬ 
torized  and  stored,  solving  for  multiple  right-hand-side 
vectors  by  back  substitution  would  be  relatively  inex¬ 
pensive.  However,  for  very  large  systems  —  such  as 
the  ones  encountered  in  CFD  —  the  matrix  dTZk/dyi 
is  never  factorized  explicitly  and  the  system  of  equa¬ 
tions  requires  an  iterative  solution  which  is  usually  as 
costly  as  solving  the  governing  equations.  When  we 
multiply  this  cost  by  the  number  of  design  variables, 
the  total  cost  for  calculating  the  sensitivity  vector  may 
become  unacceptable. 

Returning  to  the  total  sensitivity  equation  (12),  we 
observe  that  there  is  an  alternative  option  when  com¬ 
puting  the  total  sensitivity  dl /  dxn.  The  auxiliary 
vector  Tfc  can  be  obtained  by  solving  the  adjoint  equa¬ 
tions 


dl_ 
dyi  ’ 


(13) 


The  vector  d>k  is  usually  called  the  adjoint  vector  and 
is  substituted  into  equation  (12)  to  find  the  total  sensi¬ 
tivity.  In  contrast  with  the  direct  method,  the  adjoint 
vector  does  not  depend  on  the  design  variables,  xn, 
but  instead  depends  on  the  function  of  interest,  I. 

We  can  now  see  that  the  choice  of  solution  procedure 
(direct  vs.  adjoint)  to  obtain  the  total  sensitivity  (12) 
has  a  substantial  impact  on  the  cost  of  sensitivity  anal¬ 
ysis.  Although  all  the  partial  derivative  terms  are  the 
same  for  both  the  direct  and  adjoint  methods,  the 
order  of  the  operations  is  not.  Notice  that  for  any  num¬ 
ber  of  functions,  /,  we  can  compute  dyi/ dxn  once  for 
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each  design  variable  (direct  method).  Alternatively, 
for  an  arbitrary  number  of  design  variables,  we  can 
compute  d'fe  once  for  each  function  (adjoint  method). 

The  cost  involved  in  calculating  sensitivities  using 
the  adjoint  method  is  therefore  practically  indepen¬ 
dent  of  the  number  of  design  variables.  After  having 
solved  the  governing  equations,  the  adjoint  equations 
are  solved  only  once  for  each  I.  Moreover,  the  cost  of 
solution  of  the  adjoint  equations  is  similar  to  that  of 
the  solution  of  the  governing  equations  since  they  are 
of  similar  complexity  and  the  partial  derivative  terms 
are  easily  computed. 

Therefore,  if  the  number  of  design  variables  is 
greater  than  the  number  of  functions  for  which  we  seek 
sensitivity  information,  the  adjoint  method  is  compu¬ 
tationally  more  efficient.  Otherwise,  if  the  number  of 
functions  to  be  differentiated  is  greater  than  the  num¬ 
ber  of  design  variables,  the  direct  method  would  be  a 
better  choice. 

The  adjoint  method  has  been  widely  used  for  both 
structural  sensitivity  analysis  [1]  and  aerodynamic 
shape  optimization  [13,  14]. 

Aero-Structural  Adjoint  Equations 

Although  the  theory  we  have  just  presented  is  appli¬ 
cable  to  multidisciplinary  systems,  provided  that  the 
governing  equations  for  all  disciplines  are  included  in 
IZk,  we  now  explicitly  discuss  the  sensitivity  analy¬ 
sis  of  multidisciplinary  systems,  using  aero-structural 
optimization  as  an  example.  This  example  illustrates 
the  fundamental  computational  cost  issues  that  mo¬ 
tivate  our  choice  of  strategy  for  sensitivity  analysis. 
The  following  equations  and  discussion  can  easily  be 
generalized  for  cases  with  additional  disciplines. 

In  the  aero-structural  case  we  have  coupled  aerody¬ 
namic  (Ak)  and  structural  (Si)  governing  equations, 
and  two  sets  of  state  variables:  the  flow  state  vec¬ 
tor,  Wi,  and  the  vector  of  structural  displacements, 
Uj.  In  the  following  expressions,  we  split  the  vectors 
of  residuals,  states  and  adjoints  into  two  smaller  vec¬ 
tors  corresponding  to  the  aerodynamic  and  structural 
systems 


TZk' 


Ak 

Si 


Figure  4  shows  a  diagram  representing  the  coupling  in 
this  system. 


Fig.  4:  Schematic  representation  of  the  aero-structural 
system. 


Using  this  new  notation,  the  adjoint  equations  (13) 
for  the  case  of  the  aero-structural  system  can  be  writ¬ 
ten  as 


VdAk  dAk 1 

T 

r  dr  1 

dwi  duj 

4>k 

dwi 

dSi  dSi 

4>i 

— 

di 

_  dwi  duj  _ 

|_9u7J 

Note  that  the  residual  sensitivity  matrix  in  this  equa¬ 
tion  is  identical  to  that  of  the  Global  Sensitivity  Equa¬ 
tions  (GSE)  first  presented  by  Sobieski  [27].  This 
matrix,  in  addition  to  containing  the  diagonal  terms 
that  appear  when  we  solve  the  single  discipline  ad¬ 
joint  equations,  also  has  off-diagonal  terms  expressing 
the  sensitivity  of  one  discipline  to  the  state  variables  of 
the  other.  The  details  of  the  partial  derivative  terms 
in  this  matrix  and  the  right-hand  side  (for  cases  when 
I  =  Cd  and  I  =  KS)  of  this  equation  (15)  have 
been  previously  published  and  will  not  be  repeated 
here  [18,  19]. 

Since  the  factorization  of  the  full  matrix  in  the 
coupled-adjoint  equations  (15)  would  be  extremely 
costly,  our  approach  uses  an  iterative  solver,  much  like 
the  one  used  for  the  aero-structural  solution,  where  the 
adjoint  vectors  are  lagged  and  the  two  different  sets  of 
equations  are  solved  separately.  For  the  calculation 
of  the  adjoint  vector  of  one  discipline,  we  use  the  ad¬ 
joint  vector  of  the  other  discipline  from  the  previous 
iteration,  i.e.,  we  solve 


Aerodynamic  adjoint 


dSL  =_dL 

duj  1  duj 

Structural  adjoint 


dAk 

duj 


i’k, 


(17) 


where  and  cf>i  are  the  lagged  aerodynamic  and  struc¬ 
tural  adjoint  vectors  respectively.  Upon  convergence, 
the  final  result  given  by  this  system  is  the  same  as 
that  of  the  original  coupled-adjoint  equations  (15). 
We  call  this  the  lagged-coupled  adjoint  (LCA)  method 
for  computing  sensitivities  of  multidisciplinary  sys¬ 
tems.  Note  that  these  equations  look  like  the  single 
discipline  adjoint  equations  for  the  aerodynamic  and 
structural  solvers,  with  the  addition  of  forcing  terms 
on  the  right-hand  side  that  contain  the  off-diagonal 
terms  of  the  residual  sensitivity  matrix.  This  allows 
us  to  use  existing  single-discipline  adjoint  sensitivity 
analysis  methods  with  only  small  modifications.  Note 
also  that,  even  for  more  than  two  disciplines,  this  it¬ 
erative  solution  procedure  is  nothing  more  than  the 
well-known  block- Jacobi  method. 

Once  both  adjoint  vectors  have  converged,  we  can 
compute  the  final  sensitivities  of  the  objective  function 
by  using  the  following  expression 


d/  dl  dAk  dSi 
dxn  dxn+ipkdxn+<pldxn' 


(18) 


6  OF  14 


American  Institute  of  Aeronautics  and  Astronautics  Paper.  2002-5402 


which  is  the  coupled  version  of  the  total  sensitivity 
equation  (12). 

Note  again,  that  the  details  of  the  partial  derivative 
terms  in  the  LCA  equations  (16)  and  (17)  and  the 
total  sensitivity  equation  (18)  can  be  found  in  previous 
publications  [18,  19]. 

For  the  aero-structural  optimization  problem  at 
hand  the  aerodynamic  portion  is  usually  character¬ 
ized  by  a  single  objective  function  and  at  most  a  few 
aerodynamic  constraints,  but  a  large  number  of  design 
variables.  On  the  other  hand,  the  structural  portion 
of  the  optimization  problem  is  characterized  by  a  large 
number  of  constraints:  the  stress  in  each  element  of  the 
finite-element  model  cannot  exceed  the  material  yield 
stress  for  a  number  of  load  conditions.  Constrained 
gradient  optimization  methods  generally  require  that 
the  user  provide  the  gradient  of  both  the  cost  function 
and  each  nonlinear  constraint  with  respect  to  all  of 
the  design  variables  in  the  problem.  Using  the  adjoint 
approach,  the  evaluation  of  the  gradient  of  each  con¬ 
straint  would  require  an  independent  coupled  solution 
of  a  large  adjoint  system.  Since  the  number  of  struc¬ 
tural  constraints  is  similar  to  the  number  of  design 
variables  in  the  problem  (0(1O3)  or  larger),  the  useful¬ 
ness  of  the  adjoint  approach  could  be  put  in  question. 

Both  of  the  remaining  alternatives,  the  direct  and 
finite-difference  methods,  prove  overly  costly  since 
they  both  require  a  number  of  solutions  that  is  compa¬ 
rable  to  the  number  of  design  variables.  In  the  absence 
of  other  choices  that  can  efficiently  evaluate  the  gra¬ 
dient  of  a  large  number  of  constraints  with  respect  to 
a  large  number  of  design  variables,  it  is  necessary  to 
reduce  the  size  of  the  problem  either  through  a  reduc¬ 
tion  in  the  number  of  design  variables  or  through  a 
reduction  in  the  number  of  nonlinear  constraints. 

The  reason  for  the  choice  of  the  KS  functions  to 
lump  the  structural  constraints  now  becomes  clear. 
By  employing  KS  functions,  the  number  of  struc¬ 
tural  constraints  for  the  problem  can  be  reduced  from 
0(1O3)  to  just  a  few.  Since  the  KS  function  is  a  lower 
bound  envelope  of  all  the  constraint  inequalities,  this 
dramatic  reduction  in  the  number  of  structural  con¬ 
straints  can  in  theory  lead  to  conservative  designs.  In 
our  experience  and  that  of  other  researchers  [2],  the 
degree  of  conservativeness  added  by  the  KS  functions 
is  small. 


Results 

In  this  section  we  present  the  results  of  the  ap¬ 
plication  of  our  sensitivity  calculation  method  to  the 
problem  of  aero-structural  design  of  a  supersonic,  nat¬ 
ural  laminar  flow,  business  jet.  Before  describing  the 
results  of  our  design  experience,  we  present  the  aero- 
structural  analysis  framework  and  the  results  of  a  sen¬ 
sitivity  validation  study. 


Fig.  5:  Aero-structural  model  and  solution  of  the  su¬ 
personic  business  jet  configuration,  showing  a  slice  of 
the  grid  and  the  internal  structure  of  the  wing. 


Aero-Structural  Analysis 

The  coupled-adjoint  procedure  was  implemented 
as  a  module  that  was  added  to  the  aero-structural 
design  framework  previously  developed  by  the  au¬ 
thors  [25,  18].  The  framework  consists  of  an  aero¬ 
dynamic  analysis  and  design  module  (which  includes 
a  geometry  engine  and  a  mesh  perturbation  algo¬ 
rithm),  a  linear  finite-element  structural  solver,  an 
aero-structural  coupling  procedure,  and  various  pre¬ 
processing  tools  that  are  used  to  setup  aero-structural 
design  problems.  The  multi-disciplinary  nature  of  this 
solver  is  illustrated  in  Figure  5  where  we  can  see  the 
aircraft  geometry,  the  flow  mesh  and  solution,  and  the 
primary  structure  inside  the  wing. 

The  aerodynamic  analysis  and  design  module, 
SYN107-MB  [23],  is  a  multiblock  parallel  flow  solver 
for  both  the  Euler  and  the  Reynolds  Averaged  Navier- 
Stokes  equations  that  has  been  shown  to  be  accurate 
and  efficient  for  the  computation  of  the  flow  around 
full  aircraft  configurations  [26].  An  aerodynamic  ad¬ 
joint  solver  is  also  included  in  this  package  in  order 
to  perform  aerodynamic  shape  optimization  in  the  ab¬ 
sence  of  aero-structural  interaction. 

The  structural  analysis  package  is  FESMEH,  a  finite 
element  solver  developed  by  Holden  [12].  The  package 
is  a  linear  finite-element  solver  that  incorporates  two 
element  types  and  computes  the  structural  displace¬ 
ments  and  stresses  of  wing  structures.  Although  this 
solver  is  not  as  general  as  some  commercially-available 
packages,  it  is  still  representative  of  the  challenges  in¬ 
volved  in  using  large  models  with  tens  of  thousands  of 
degrees  of  freedom.  High-fidelity  coupling  between  the 
aerodynamic  and  the  structural  analysis  programs  is 
achieved  using  a  linearly  consistent  and  conservative 
scheme  [25,  7]. 

The  structural  model  of  the  wing  can  be  seen  in 
Figure  5  and  is  constructed  using  a  wing  box  with  six 
spars  evenly  distributed  from  15%  to  80%  of  the  chord. 
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Ribs  are  distributed  along  the  span  at  every  tenth  of 
the  semispan.  A  total  of  640  finite  elements  were  used 
in  the  construction  of  this  model.  Appropriate  thick¬ 
nesses  of  the  spar  caps,  shear  webs,  and  skins  were 
chosen  according  to  the  expected  loads  for  this  design. 

Aero-Structural  Sensitivity  Validation 


Design  variable,  n 

Fig.  6:  Sensitivities  of  the  drag  coefficient  with  respect 
to  shape  perturbations. 

To  gain  confidence  in  the  effectiveness  of  the  aero- 
structural  coupled-adjoint  sensitivities  for  use  in  de¬ 
sign  optimization,  we  must  ensure  that  the  values  of 
the  gradients  are  accurate.  For  validation  purposes, 
we  use  four  sets  of  sensitivities.  Results  from  the  ad¬ 
joint  method  are  compared  to  the  exact  discrete  value 
of  these  sensitivities  using  the  complex-step  derivative 
approximation  [20,  21]. 

In  this  sensitivity  study  two  different  functions  are 
considered:  the  aircraft  drag  coefficient,  Cd,  and  the 
KS  function  (6) .  The  sensitivities  of  these  two  quanti¬ 
ties  with  respect  to  both  OML  shape  design  variables 
and  structural  design  variables  are  computed  and  dis¬ 
cussed. 

Cd  with  Respect  to  OML  Variables 

The  values  of  the  aero-structural  sensitivities  of  the 
drag  coefficient  with  respect  to  shape  perturbations 
are  shown  in  Figure  6.  The  ten  shape  perturbations 
were  chosen  to  be  Hicks-Henne  bumps  distributed 
chordwise  on  the  upper  surface  of  two  adjacent  air¬ 
foils  around  the  quarter  span.  The  plot  shows  very 
good  agreement  between  the  coupled-adjoint  and  the 
complex-step  results,  with  an  average  relative  error 
between  the  two  of  only  3.5%.  Note  that  all  these  sen¬ 
sitivities  are  total  sensitivities  in  the  sense  that  they 
account  for  the  coupling  between  aerodynamics  and 
structures. 

To  verify  the  need  for  taking  the  coupling  into  ac¬ 
count,  the  same  set  of  sensitivities  was  calculated  for 
fixed  structural  displacements,  where  the  displacement 


field  is  frozen  after  the  aero-structural  solution.  This 
is  similar  to  assuming  that  the  wing,  after  the  ini¬ 
tial  aeroelastic  deformation,  is  held  rigid  as  far  as  the 
computation  of  sensitivities  is  concerned.  The  cal¬ 
culation  of  the  sensitivities  only  takes  into  account 
variations  related  to  the  aerodynamics.  Figure  6  shows 
that  the  single-system  sensitivities  exhibit  significantly 
lower  magnitudes  and  even  opposite  signs  for  many  of 
the  design  variables,  when  compared  with  the  coupled 
sensitivities.  The  use  of  single-discipline  sensitivities 
would  clearly  lead  to  erroneous  design  decisions. 


Design  variable,  n 

Fig.  7:  Sensitivities  of  the  drag  coefficient  with  respect 
to  structural  thicknesses. 

Cd  with  Respect  to  Thickness  Variables 

Figure  7  also  shows  the  sensitivity  of  the  drag  co¬ 
efficient,  this  time  with  respect  to  the  thicknesses  of 
five  skin  groups  and  five  spar  groups  distributed  along 
the  span.  The  agreement  in  this  case  is  even  better; 
the  average  relative  error  is  only  1.6%.  Even  though 
these  are  sensitivities  with  respect  to  internal  struc¬ 
tural  variables  that  do  not  modify  the  jig  OML,  the 
non-zero  values  in  Figure  7  demonstrate  that  coupled 
sensitivity  analysis  is  needed. 

KS  with  Respect  to  OML  and  Thickness  Variables 

The  sensitivities  of  the  KS  function  with  respect  to 
the  two  sets  of  design  variables  described  above  are 
shown  in  Figures  8  and  9.  The  results  show  that 
the  coupled-adjoint  sensitivities  are  extremely  accu¬ 
rate,  with  average  relative  errors  of  2.9%  and  1.6%. 
In  Figure  9  we  observe  that  the  sensitivity  of  the  KS 
function  with  respect  to  the  first  structural  thickness 
is  much  higher  than  the  remaining  sensitivities.  This 
markedly  different  magnitude  is  due  to  the  fact  that 
this  particular  structural  design  variable  corresponds 
to  the  thickness  of  the  top  and  bottom  skins  of  the 
wing  bay  closest  to  the  root,  where  the  stress  is  the 
highest. 

The  sensitivities  of  the  KS  function  for  fixed  loads 
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Fig.  8:  Sensitivities  of  the  KS  function  with  respect  to 
shape  perturbations. 


Design  variable,  n 

Fig.  9:  Sensitivities  of  the  KS  function  with  respect  to 
structural  thicknesses. 

are  also  shown  in  Figures  8  and  9.  Using  the  complex- 
step  method,  these  sensitivities  were  calculated  by 
calling  only  the  structural  solver  after  the  initial  aero- 
structural  solution.  The  approach  is  equivalent  to  us¬ 
ing  just  equations  (17,  18)  without  the  partial  deriva¬ 
tives  of  Ak ■  The  difference  in  these  sensitivities  when 
compared  to  the  coupled  ones  is  not  as  dramatic  as  in 
the  fixed  displacements  case  shown  in  Figure  6,  but  it 
is  still  significant. 

Computational  Efficiency 

The  cost  of  calculating  a  gradient  vector  using  ei¬ 
ther  the  finite-difference  or  the  complex-step  methods 
is  expected  to  be  linearly  dependent  on  the  number  of 
design  variables.  This  expectation  is  confirmed  in  Fig¬ 
ure  10  where  the  gradient  calculation  times  are  shown 
for  an  increasing  number  of  design  variables.  The  time 


axis  is  normalized  with  respect  to  the  time  required 
for  a  single  aero-structural  solution  (98  seconds  on  9 
processors  of  an  SGI  Origin  2000).  The  symbols  indi¬ 
cate  timings  measured  from  actual  calculations:  2, 000 
design  variables  were  tested  for  the  coupled-adjoint 
method. 
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Fig.  10:  Computational  time  vs.  number  of  design  vari¬ 
ables  for  finite  differencing,  complex  step  and  coupled 
adjoint.  The  time  is  normalized  with  respect  to  the 
time  required  for  one  aero-structural  solution. 

The  cost  of  a  finite-difference  gradient  evaluation 
can  be  linearly  approximated  by  the  equation  1.0  + 
0.38  x  Nx,  where  Nx  is  the  number  of  design  variables. 
Notice  that  one  might  expect  this  method  to  incur  a 
computational  cost  equivalent  to  one  aero-structural 
solution  per  additional  design  variable.  The  cost  per 
additional  design  variable  is  lower  than  unity  because 
each  additional  aero-structural  calculation  does  not 
start  from  a  uniform  flow-field  initial  condition,  but 
from  the  previously  converged  solution. 

The  same  applies  to  the  cost  of  the  complex-step 
method.  Because  the  function  evaluations  require 
complex  arithmetic,  the  cost  of  the  complex  step 
method  is,  on  average,  2.4  times  higher  than  that  of 
finite  differencing.  However,  this  cost  penalty  is  worth¬ 
while  since  there  is  no  need  to  find  an  acceptable  step 
size  a  priori ,  as  is  the  case  for  finite-difference  approx¬ 
imations  [20,  21]. 

The  cost  of  computing  sensitivities  using  the 
coupled-adjoint  procedure  is  in  theory  independent  of 
the  number  of  variables.  Using  our  implementation, 
however,  some  of  the  partial  derivatives  in  the  total 
sensitivity  equation  (18)  are  calculated  using  finite  dif¬ 
ferences  involving  mesh  perturbations  and  therefore, 
there  is  a  small  dependence  on  the  number  of  vari¬ 
ables.  The  line  representing  the  cost  of  the  coupled 
adjoint  in  Figure  10  has  a  slope  of  0.01  which  is  be¬ 
tween  one  and  two  orders  of  magnitude  less  than  the 
slope  for  the  other  two  lines. 

The  constant  terms  for  the  straight  lines  in  Figure  10 
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represent  the  upfront  cost  of  each  procedure  when  no 
sensitivities  are  required.  For  the  finite-difference  case, 
this  is  equivalent  to  one  aero-structural  solution,  and 
hence  the  constant  is  1.0.  When  performing  the  aero- 
structural  solution  using  complex  arithmetic,  the  cost 
rises  to  2.1  times  the  real  arithmetic  solution.  The 
cost  of  computing  the  coupled-adjoint  vectors  (with¬ 
out  computing  the  gradients)  is  3.4.  This  cost  includes 
the  aero-structural  solution,  which  is  necessary  before 
solving  the  adjoint  equations,  and  hence  the  aero- 
structural  adjoint  computation  alone  incurs  a  cost  of 
2.4.  Since  the  solution  of  the  aero-structural  adjoint 
equation  incurs  a  computational  cost  similar  to  that  of 
the  aero-structural  analysis,  the  equivalent  of  1.4  aero- 
structural  analyses  is  spent  in  the  computation  of  the 
forcing  terms  introduced  in  equations  (16)  and  (17). 
In  particular,  the  last  term  in  equation  (17)  requires 
that  the  full  CFD  mesh  be  perturbed  for  each  of  the 
degrees  of  freedom  on  the  surface  of  the  CSM  mesh. 
The  calculation  of  this  term  can  increase  the  compu¬ 
tational  cost  for  very  large  structural  models. 

The  main  conclusion  still  remains  that  the  cost  of 
computing  sensitivities  with  respect  to  hundreds  or 
even  thousands  of  variables  is  acceptable  when  using 
the  coupled-adjoint  approach,  while  it  is  impractical 
to  use  finite-differences,  the  complex-step  method,  or 
the  direct  method  for  such  a  large  number  of  design 
variables. 

Aero-Structural  Design 

The  objective  in  this  optimization  is  to  solve  the  op¬ 
timization  problem  that  we  previously  described,  i.e., 


minimize 

1  = 

z  o.C jj  +  j 

xn  G  M™ 
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In  our  example  the  value  of  Cd  corresponds  to  that  of 
the  cruise  condition,  which  has  a  target  lift  coefficient 
of  0.1.  The  structural  stresses,  in  the  form  of  the  KS 
function,  correspond  to  a  single  maneuver  condition, 
for  which  Clt  =0.2. 

All  optimization  work  is  carried  out  using  the  non¬ 
linear  constrained  optimizer  NPSOL  [10].  Euler  calcu¬ 
lations  are  performed  on  a  wing-body  36-block  mesh 
that  is  constructed  from  the  decomposition  of  a  193  x 
33  x  49  C-H  mesh.  During  the  process  of  the  optimiza¬ 
tion,  all  flow  evaluations  are  converged  to  5.3  orders  of 
magnitude  of  the  average  density  residual  and  the  Cl 
constraint  is  achieved  to  within  1CD6. 

In  order  to  parameterize  the  shape  of  the  aircraft,  we 
have  chosen  sets  of  design  variables  that  apply  to  both 
the  wing  and  the  fuselage.  The  wing  shape  is  modified 
by  the  design  optimization  procedure  at  six  defining 
stations  uniformly  distributed  from  the  side-of-body 


to  the  tip  of  the  wing.  The  shape  modifications  of 
these  defining  stations  are  linearly  lofted  to  a  zero 
value  at  the  previous  and  next  defining  stations.  On 
each  defining  station,  the  twist,  the  leading  and  trail¬ 
ing  edge  camber  distributions,  and  five  Hicks-Henne 
bump  functions  on  both  the  upper  and  lower  surfaces 
are  allowed  to  vary.  The  leading  and  trailing  edge  cam¬ 
ber  modifications  are  not  applied  at  the  first  defining 
station.  This  yields  a  total  of  76  OML  design  variables 
on  the  wing.  Planform  modifications,  which  are  per¬ 
mitted  by  our  software,  were  not  used  in  the  present 
calculations.  Planform  optimization  is  only  meaning¬ 
ful  if  additional  disciplines  and  constraints  are  taken 
into  account. 

The  shape  of  the  fuselage  is  parameterized  in  such  a 
way  that  its  camber  is  allowed  to  vary  while  the  total 
volume  remains  constant.  This  is  accomplished  with 
9  bump  functions  evenly  distributed  in  the  streamwise 
direction  starting  at  the  10%  fuselage  station.  Fuselage 
nose  and  trailing  edge  camber  functions  are  added  to 
the  fuselage  camber  distribution  in  a  similar  way  to 
what  was  done  with  the  wing  sections. 

The  structural  sizing  is  accomplished  with  10  design 
variables,  which  correspond  to  the  skin  thicknesses  of 
the  top  and  bottom  surfaces  of  the  wing.  Each  group 
is  formed  by  the  plate  elements  located  between  two 
adjacent  ribs.  All  structural  design  variables  are  con¬ 
strained  to  exceed  a  pre-specified  minimum  gage  value. 

The  complete  configuration  is  therefore  parameter¬ 
ized  with  a  total  of  97  design  variables.  As  mentioned 
in  an  earlier  section,  the  cost  of  aero-structural  gra¬ 
dient  information  using  our  coupled-adjoint  method  is 
effectively  independent  of  the  number  of  design  vari¬ 
ables:  in  more  realistic  full  configuration  test  cases 
that  we  are  about  to  tackle,  500  or  more  design  vari¬ 
ables  will  be  necessary  to  describe  the  shape  variations 
of  the  configuration  (including  nacelles,  diverters,  and 
tail  surfaces)  and  the  sizing  of  the  structure. 

The  initial  application  of  our  design  methodology  to 
the  aero-structural  design  of  a  supersonic  business  jet 
is  simply  a  proof-of-concept  problem  meant  to  validate 
the  sensitivities  obtained  with  our  method.  Current 
work  is  addressing  the  use  of  multiple  realistic  load 
conditions,  dynamic  loads,  aeroelastic  constraints,  and 
the  addition  of  diverters,  nacelles,  and  empennage. 

In  the  present  design  case,  we  use  a  =  104,  /3  = 
3.226  x  10-3.  Note  that  the  scalars  that  multiply  the 
structural  weight,  W,  and  the  coefficient  of  drag,  Cd, 
reflect  the  correct  trade-off  between  drag  and  weight 
that  was  previously  mentioned,  i.e.  that  one  count  of 
drag  is  worth  310  pounds  of  weight. 

Figure  11  shows  the  evolution  of  this  aero-structural 
design  case  for  successive  major  design  iterations.  The 
figure  shows  the  values  of  the  coefficient  of  drag  (in 
counts),  the  wing  structural  weight  (in  lbs),  and  the 
value  of  the  KS  function.  Note  that  the  structural  con¬ 
straints  are  satisfied  when  the  KS  function  is  positive. 
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Fig.  11:  Convergence  history  of  the  aero-structural  optimization. 
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Because  of  the  approximate  nature  of  the  KS  function, 
all  structural  constraints  may  actually  be  satisfied  for 
small  but  negative  values  of  the  KS  function. 

The  baseline  design  is  feasible,  with  a  cruise  drag 
coefficient  of  74.04  counts  and  a  structural  weight  of 
9,  285  lbs.  The  KS  function  is  slightly  positive  indi¬ 
cating  that  all  stress  constraints  are  satisfied  at  the 
maneuver  condition.  In  the  first  two  design  iterations, 
the  optimizer  takes  large  steps  in  the  design  space, 
resulting  in  a  drastic  reduction  in  both  Cd  and  W. 
However,  this  also  results  in  a  highly  infeasible  de¬ 
sign  that  exhibits  maximum  stresses  2.1  times  the  yield 
stress  of  the  material.  After  these  initial  large  steps, 
the  optimizer  manages  to  decrease  the  norm  of  the 
constraint  violation.  This  seems  to  have  been  accom¬ 
plished  by  increasing  the  structural  skin  thicknesses, 
since  the  weight  increases  while  the  drag  is  further  re¬ 
duced.  Towards  major  iteration  10,  there  is  no  visible 
progress  for  several  iterations  while  the  design  remains 
infeasible.  A  large  step  is  taken  in  iteration  13  that 
results  in  a  sudden  increase  in  feasibility  accompanied 
by  an  equally  sudden  increase  in  Cd-  The  optimizer 
has  established  that  the  only  way  to  obtain  a  feasible 
design  is  by  increasing  the  wing  thickness  (with  the 
consequent  increases  in  Cd  and  weight)  and  the  struc¬ 
tural  thicknesses.  From  that  point  on,  the  optimizer 
rapidly  converges  to  the  optimum.  After  43  major  iter¬ 
ations,  the  KS  constraint  is  reduced  to  0(1O~4)  and  all 
stress  constraints  are  satisfied.  The  aero-structurally 
optimized  result  has  Cd  =  0.006922  and  a  total  wing 
structure  weight  of  5, 546  lbs. 

Visualizations  of  the  baseline  and  optimized  config¬ 
urations  are  shown  in  Figures  12  and  13.  Measures 
of  performance  and  feasibility  are  written  in  the  first 
section  of  Table  1.  The  left  halves  of  Figures  12  and  13 


show  the  surface  flow  density  distribution  with  the 
corresponding  structural  deflections  at  the  cruise  con¬ 
dition  for  both  the  initial  and  optimized  designs.  The 
right  halves  show  an  exploded  view  of  the  stress  dis¬ 
tribution  on  the  structure  (spar  caps,  spar  shear  webs, 
and  skins,  from  top  to  bottom)  at  the  Cl  =  0.2  maneu¬ 
ver  condition.  From  these  Figures  one  can  appreciate 
that  not  only  have  the  surface  density  distributions 
changed  substantially  at  the  cruise  point,  but  so  have 
the  element  stresses  at  the  maneuver  condition.  In 
fact,  as  expected  from  a  design  case  with  a  single 
load  condition,  the  optimized  structure  is  nearly  fully- 
stressed,  except  in  the  outboard  sections  of  the  wing, 
where  the  minimum  gage  constraints  are  active.  It  is 
also  worth  noting  that  about  half  of  the  improvement 
in  the  Cd  of  the  optimized  configuration  results  from 
drastic  changes  in  the  fuselage  shape:  both  front  and 
aft  camber  have  been  added  to  distribute  the  lift  more 
evenly  in  the  streamwise  direction  in  order  to  reduce 
the  total  lift-dependent  wave  drag. 

A  total  of  50  major  design  iterations  including  aero- 
structural  analyses,  coupled  adjoint  solutions,  gradient 
computations,  and  line  searches  were  performed  in 
approximately  20  hours  of  wall  clock  time  using  18 
processors  of  an  SGI  Origin  3000  system  (R12000, 
400  Mhz  processors).  Since  these  are  not  the  fastest 
processors  currently  available  we  feel  confident  that 
much  larger  models  can  be  optimized  with  overnight 
turnaround  in  the  near  future. 

Comparison  With  Sequential  Optimization 

The  usefulness  of  a  coupled  aero-structural  opti¬ 
mization  method  can  only  be  measured  in  compari¬ 
son  with  the  results  that  can  be  obtained  using  cur¬ 
rent  state-of-the-art  practices.  In  the  case  of  aero- 
structural  design,  the  typical  approach  is  to  carry  out 
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aerodynamic  shape  optimization  with  artificial  thick¬ 
ness  constraints  meant  to  represent  the  effect  of  the 
structure,  followed  by  structural  optimization  with  a 
fixed  OML.  It  is  well  known  that  sequential  optimiza¬ 
tion  cannot  be  guaranteed  to  convergence  to  the  true 
optimum  of  a  coupled  system.  In  order  to  determine 
the  difference  between  the  optima  achieved  by  fully- 
coupled  and  sequential  optimizations,  we  have  also 
carried  out  one  cycle  of  sequential  optimization  within 
our  analysis  and  design  framework. 

To  prevent  the  optimizer  from  thinning  the  wing  to 
an  unreasonable  degree  during  the  aerodynamic  shape 
optimization,  5  thickness  constraints  are  added  to  each 
of  the  6  defining  stations  for  a  total  of  30  linear  con¬ 
straints.  These  constraints  are  such  that,  at  the  points 
where  they  are  applied,  the  wing  box  is  not  allowed  to 
get  any  thinner  than  the  original  design. 

After  the  process  of  aerodynamic  shape  optimiza¬ 
tion  is  completed,  the  initial  Cd  has  decreased  to 
Cd  =  0.006992  as  can  be  seen  in  the  second  half  of  Ta¬ 
ble  1.  After  fixing  the  OML,  structural  optimization 
is  performed  using  the  maneuver  loads  for  the  baseline 
configuration  at  Cl  =0.2.  The  optimized  structural 
design  reduces  the  empty  weight  of  the  structure  to 
6,  567  lbs. 

We  can  now  compare  the  results  of  the  fully  cou¬ 
pled  optimization  in  the  previous  section  and  the  out¬ 
come  of  the  process  of  sequential  optimization.  The 
differences  are  clear:  the  coupled  aero-structural  opti¬ 
mization  was  able  to  achieve  an  optimized  design  with 
a  weight  of  only  5,  546  lbs  when  compared  with  the 
6,  567  lbs  of  the  sequential  optimization.  This  repre¬ 
sents  a  relative  difference  of  nearly  16%.  Although  it 
is  barely  distinguishable,  the  coefficient  of  drag  for  the 
coupled  optimization  approach  is  slightly  lower  than 
for  sequential  optimization. 

Finally,  notice  that  since  sequential  optimization  ne¬ 
glects  the  aero-structural  coupling  in  the  computation 
of  maneuver  loads,  there  is  no  guarantee  that  the 
resulting  design  will  be  feasible.  In  fact,  the  aero- 
structural  analysis  shows  that  the  value  of  the  KS 
function  is  slightly  negative. 

Conclusions 

A  methodology  for  coupled  sensitivity  analysis  of 
high-fidelity  aero-structural  systems  was  presented. 
The  sensitivities  computed  by  the  lagged-coupled- 
adjoint  method  were  compared  to  sensitivities  given  by 
the  complex-step  derivative  approximation  and  shown 
to  be  extremely  accurate,  having  an  average  relative 
error  of  2%.  Moreover,  significant  differences  in  the 
values  and  signs  of  the  sensitivities  were  found  when 
aero-structural  values  were  compared  to  rigid  ones. 

In  realistic  aero-structural  design  problems  with 
hundreds  of  design  variables,  there  is  a  considerable  re¬ 
duction  in  computational  cost  when  using  the  coupled- 
adjoint  method  as  opposed  to  either  finite-differences 


or  the  complex-step  approaches.  This  improvement  is 
due  to  the  fact  that  the  cost  associated  with  the  ad¬ 
joint  method  is  practically  independent  of  the  number 
of  design  variables. 

Sensitivities  computed  using  the  presented  method¬ 
ology  were  successfully  used  to  optimize  the  design  of  a 
supersonic  business  jet  that  was  parameterized  with  a 
large  number  of  aerodynamic  and  structural  variables. 
The  outcome  of  this  optimization  was  compared  with 
the  traditional  method  of  sequential  optimization  and 
it  was  found  to  improve  the  structural  weight  by  an 
additional  16%. 
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Cd  (counts) 

KS 

^"max/  ^yield 

Weight  (lbs) 

Objective 

Integrated  approach 

Baseline 

73.95 

1.15  x  lO^1 

0.87 

9,285 

Optimized 

69.22 

-2.68  x  10"4 

0.98 

5,546 

87.12 

Sequential  approach 

Aerodynamic  optimization 
Baseline 

74.04 

Optimized 

Structural  optimization 

69.92 

Baseline 

1.02  x  10"1 

0.89 

9,285 

Optimized 

1.45  x  10"8 

0.98 

6,567 

Aero-structural  analysis 

69.95 

-9.01  x  10"3 

0.99 

6,567 

91.14 

Table  1:  Comparison  between  the  integrated  and  sequential  approaches  to  aero-structural  optimization. 
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Derail  y 


Fig.  12:  Baseline  configuration  for  the  supersonic  business  jet  showing  surface  densities  at  the  cruise  condition 
and  structural  stresses  at  the  maneuver  condition.  The  density  is  normalized  by  the  freestream  value  and  the 
von  Mises  stresses  are  normalized  by  the  material  yield  stress. 
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Fig.  13:  Optimized  configuration  for  the  supersonic  business  jet. 
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